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We show that the effect of stress on the reproductive fitness of noisy cell populations can be 
modelled as first-passage time problem, and demonstrate that even relatively short-lived fluctuations 
in gene expression can ensure long-term survival of a drug-resistant population. We examine how this 
effect contributes to the development of drug-resistant cancer cells, and demonstrate that permanent 
immunity can arise independently of mutations. 



Gene expression is a stochastic process that enables 
genetically identical cells in the same environment to 
exhibit phenotypic variation [1-3 . This noise-induced 
non-genetic (epigenetic) variability can be beneficial to 
cell populations experiencing acute stress by providing a 
temporary basis for natural selection [3H7]. 

Experimental observations suggest that gene expres- 
sion is inherently associated with 'epigenetic memory', 
defined by the fluctuation relaxation time of a gene prod- 
uct within a cell lineage. In human lung cancer cells, this 
relaxation time can be as long as four generations [8 . 

Brock et al. [9] recently argued that epigenetic mem- 
ory might accelerate tumour progression by contribut- 
ing to the development of drug-resistant cancer cells. In 
this hypothesis, phenotypic variability from the noisy ex- 
pression of gene X that confers resistance renders some 
cells (and their offspring) temporarily insensitive to the 
drug, thereby increasing the probability of acquiring a 
mutation conferring permanent immunity. In the present 
work, we develop a minimal model to study this phe- 
nomenon quantitatively. 

To study how gene expression noise impacts the dy- 
namics of isogenic cell populations under stress, we define 
the reproductive fitness (W) as the number of offspring 
produced in the presence of the stressor (i.e. a drug) rel- 
ative to that produced in its absence. For simplicity, we 
assume that all cells produce offspring at the same rate 
in the absence of the drug, and define the generation 
time (tp) as the time it takes for each cell to reproduce 
once. We set the generation time as unit time and report 
all time-scales relative to tp. We also assume that cells 
carry the gene X conferring drug resistance when its ex- 
pression level x is sufficiently high, and that this gene is 
expressed stochastically in individual cells. 

The effects of gene expression noise on populations un- 



der stress have previously been analyzed to explain why 
certain genes have high expression noise [5H7]. In these 
analyses, the dependency between gene expression and 
reproductive fitness was defined by the integral 

W(t) = J w(x)p x (x,t)dx J (1) 

where p x (x,t) is the probability distribution function 
(PDF) describing the concentration (x) of the gene prod- 
uct across the population, and w(x) is the microscopic 
fitness function describing the effect of the drug on the 
fitness of cells with a given expression level. The basic 
concept is illustrated in Fig. [TJa) using a model where 
w(x) is described by the Heaviside step function, such 
that cells are unable to reproduce if their expression level 
is below a critical value, w(x < x c ) = 0, and unaffected 
by the drug otherwise, w(x > x c ) = 1. In this case, 
previous theoretical work [5 - 7 concluded that high gene 
expression noise is beneficial at high drug-doses, since the 
fraction of cells expressing above a reproductive thresh- 
old x c increases with the width of the initial expression 
distribution [Fig. [TJa)]. However, because p x (x, i) is as- 
sumed fixed at the time of drug treatment, this conclu- 
sion is valid only for instantaneous selection effects. The 
analysis of prolonged stress exposure necessitates an ap- 
proach where selection, inheritance and gene expression 
dynamics all contribute to the evolution of the popula- 
tion. 

Population survival during prolonged drug exposure is 
a first-passage time problem. In the absence of muta- 
tions conferring permanent immunity, cells that survive 
the initial selection will eventually succumb to the drug 
since they cannot maintain high expression indefinitely. 
Consider a sub-population of cells with the same level 
of x above x c [Fig. [lib)]. The time interval in which 
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a given cell can reproduce is the first-passage (or so- 
journ) time ts(x), where the threshold x c represents an 
absorbing barrier. Although cells are initially identical, 
the expression of the drug-resistance gene evolves differ- 
ently in different cells, and the time to reach the repro- 
ductive threshold is a random variable described by the 
first-passage time distribution ps(%,ts) (Fig. [TJb), In- 
set). Since only cells with ts(x) > to reproduce, w(x) in 
Eq. [I] is given by 
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and the overall fitness of the population at time t can be 
written as 



W(t) 
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p s (x, t s )dt' s ) p x (x, t)dx. (3) 



The population fitness in Equation ([3| has an explicit 
solution only in special cases. Previous analyses [5-7 
circumvented this problem, in part, by focusing on initial 
selection effects (t —> 0). However, even in this limit, it is 
also necessary to assume that all cells above the threshold 
contribute to fitness (i.e., w{x) = 1 for x > x c ). 

To investigate more general cases, we used the 
Ornstein-Uhlenbeck (OU) process to model the level of 



gene expression in individual cells [10 j. This process can 
be described by the Langevin equation 



dx(t) 
dt 
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where c and r are the diffusion constant and the relax- 
ation time, respectively, and £ t is Gaussian white noise 
((&) = 0, (66') = S(t - t')) QI]. The steady-state PDF 
of the OU process is a Gaussian distribution with mean 
\i and variance a 2 = cr/2. Without loss of generality, we 
set fi = and use the fluctuation time-scale r to model 
the time-scale of epigenetic memory. 

The fluctuation time-scale of gene expression has been 
determined experimentally in human lung cancer cells 
in terms of the 'mixing time' r m , defined as the lag 
where the autocorrelation function has decreased by 
50% |B]. The mixing time for the stationary OU pro- 
cess is r m = rm(2). The measured values of r m varied 
between 0.5 to 3.0 generations for different genes, corre- 
sponding to values of r between 0.7 to 4.0 generations 
for the OU process. 

First, we examined the effect of drug treatment on re- 
productive fitness after one generation time when the ab- 
sorbing barrier is located at x c = 0. In this case, the 
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FIG. 1. Epigenetic effects on a cell population exposed to 
stress, (a) Schematic of instantaneous selection effects, (b) 
Schematic of generalized model, (c) Reproductive fitness at 
the time of first division W(t = to) after the application of a 
stress (at x c = 0) as a function of r or a 2 for fixed a 2 or r, 
respectively. Analytical curves (solid lines) were obtained via 
numerical solution of Eq. (J3|. (d) W(t = as a function of 
x c for high and low a 2 , r and a 2 are scaled by £d- Dashed 
lines represent results obtained from Eq. or equivalently 
Eq. (|3| in the limit r — » oo. 




FIG. 2. Effect of epigenetic memory r on drug resistance at 
various timescales. (a) Top and bottom plots show population 
distributions corresponding respectively to short (r = 0.5) 
and long (r = 10) epigenetic memory and show the fraction of 
drug resistant cells (i.e., cells with x > x c ) after acute (t = 0) 
and a prolonged (t = 10) drug exposures (single realization of 
10 5 cells), (b) W as a function of t for various values of r. t, 
t and a 2 are scaled by to- 
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first-passage time PDF for x > x c is given by 



Ps(x,t s ) 



: exp 
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We evaluated the effects of varying the time-scale of epi- 
genetic memory and the noise amplitude by numerical 
integration of Eq. (J3|, using the steady-state OU distri- 
bution to describe the initial gene expression distribution. 
Figure [IJc) shows the results for fixed noise (a 2 = 1) and 
variable r, and fixed time-scale (r = 2) and variable a 2 . 

The time-scale of epigenetic memory significantly af- 
fects 'acute' reproductive fitness, even for very long fluc- 
tuation relaxation times. For example, when r = 20, 
W is reduced to 0.4, compared with the value of 0.5 ob- 
tained (irrespectively of the noise amplitude) in the per- 
manent epigenetic memory limit r — »> oo [Fig. [TJc)]. For 
r = 2, the reproductive fitness is approximately 0.2, and 
the majority of cells starting with x > x c are unable to 
maintain above-threshold gene expression long enough to 
reproduce. In this case, the acute reproductive fitness re- 
mains constant, presumably because changing the noise 
amplitude for x c = does not change the fraction of cells 
with x > x c . 

To examine cases where x c > 0, it is necessary to use 
numerical simulations since a general closed-form solu- 
tion of the first-passage time PDF is not available. For 
this purpose, we employed a population simulation algo- 
rithm [13] in which gene expression in each of N individ- 
ual cells, Xi(t) for i = 1,...,7V, is obtained by solving 
Eq. [4] numerically [14 . In these simulations (20 realiza- 
tions of 10 4 cells unless indicated otherwise), cell divi- 
sion occurs when a deterministic cell cycle 'clock', which 
is reset at each division, reaches tp. Each cell keeps 
track of the time since its birth and can only advance its 
clock if they maintain gene expression above the thresh- 
old. Moreover, cells where xi(t) < Xq are assumed to 
be fixed and unable to change their expression level (i.e., 
r = 00). Simulations were initiated by assigning, to each 
cell, random initial values of gene expression and the cell 
cycle clock from the steady-state distribution of the OU 
process and a uniform distribution [0 : to], respectively. 

Numerical calculations of fitness for x c > identified r 
as a critical determinant of population survival. Specifi- 
cally, the fitness of a population with low gene expression 
noise can be greater than that of a population with high 
noise if the fluctuation relaxation time is sufficiently long. 
We observed this in simulations, shown in Fig.JTJeL), with 
an increased threshold x c for fixed time-scales (r = 2 or 
r = 5) and two different fluctuation amplitudes (a 2 = 1 
or a 2 = 10). When the two populations had the same 
finite value of r, we observed that increased gene ex- 
pression noise always provides a fitness benefit (data not 
shown). However, as expected from Fig.[TJc), incorporat- 



ing stochastic gene expression dynamics (i.e., finite val- 
ues of r) generally yields a significant reduction in fitness 
compared to the asymptotic permanent memory limit. 
The magnitude of this reduction is sensitive to both the 
value of the threshold x c and the value of r. This is il- 
lustrated in Fig. [ljd) where the fitness of the high noise 
population is greater than the low noise population only 
when the value of x c is sufficiently high. 

In our second case, we analyzed the long-term effects 
of varying the time-scale of epigenetic memory on popu- 
lation dynamics and reproductive fitness. For simplicity, 
we focus on the case where x c = and noise is fixed 
(a 2 = 1). Figure [2^a) shows representative gene expres- 
sion distributions obtained after 10 generation times for 
short- and long-term epigenetic memory. When the fluc- 
tuation time-scale is short (r = 0.5, top panel), the num- 
ber of cells that may reproduce (i.e., cells with Xi(t) > x c ) 
is reduced over time since, on average, cells reach the ab- 
sorbing barrier faster than they reproduce. Correspond- 
ingly, given enough time, the population will go extinct. 
This is not the case when memory is long (r = 10, bottom 
panel) and the birth rate exceeds the rate of loss at the 
absorbing barrier. In addition, the mode of gene expres- 
sion distribution shifts to higher values, in resemblance 
of experimental observations [7]. 

Relatively short-term epigenetic memory can result in 
permanent drug resistance even in the absence of mu- 
tations. This is illustrated in Fig. [2^b) , which shows 
how the reproductive fitness of populations with differ- 
ent memory time-scales evolve over time. In populations 
with long-term memory (e.g., r = 5, 10, or 00), the num- 
ber of cells that may reproduce increases steady over time 
and settles in a steady-state where more than half of 
them reproduce every generation time (i.e., W(t) > 0.5). 
Importantly, populations with memory at intermediate 
time-scales (e.g., r = 1.5,2, or 3) may retain long-term 
viability and finite rates of reproductive fitness. Because 
the simulations involve finite populations, the outcome of 
a given realization cannot always be predicted. For ex- 
ample, when r = 1.5, a viable population was observed 
to develop in 29% of the simulations while the population 
went extinct in the remaining 71% of simulations. While 
populations with short memory (e.g., r = 0.5 or 1) even- 
tually go extinct, several cell cycles were needed for the 
drug to fully effect all cells. 

In the third and final case, we investigated the added 
effect of genetic mutations on the development of drug- 
resistance. A central element of the Brock et al. hypoth- 
esis is that temporary drug resistance due to slow fluc- 
tuations in gene expression may contribute to tumour 
development by increasing the overall probability that 
some cells acquire a mutation conferring permanent im- 
munity. To model this scenario, we allowed each cell with 
an expression level above x c the chance to mutate once 
per generation time. We denote this probability Pm- If 
a cell acquired the mutation, it and its offspring were 
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permanently resistant to the drug, and the survival of a 
continuously growing population inevitable. 

We first investigated the added effect of mutations on 
the re-emergence of a cancerous tumour under constant 
drug treatment. In these simulations, we chose x c such 
that the drug instantaneously removed 95% of the popu- 
lation, and measured the time it took for the remaining 
cells to double in number. Figure^ a) shows the depen- 
dency of this doubling time on r when Pm is equal to 
0.01 and 0.1. These mutation rates are unrealistically 
high biologically and were chosen to illustrate the effect 
of epigenetic memory in an extreme limit. 

As expected, increasing the mutation probability sig- 
nificantly reduces the doubling time when the gene ex- 
pression fluctuations are short-lived. Unexpected, how- 
ever, the value of r beyond which mutations do not have 
an additional effect is remarkably short despite the unre- 
alistically high mutation rates. Specifically, the doubling 
time is more or less unaffected by Pm when r is roughly 
above 4 generations, corresponding to the upper range of 
mixing times observed experimentally [8 . 

We confirmed our results using a semi-realistic model 
of gene expression noise [15] where proteins are synthe- 
sized in irregular bursts at irregular intervals (Fig. |3ja), 
Inset). We also tested the effect of replacing the fitness 
threshold with a more realistic sigmoidal fitness function 
and found no qualitative difference (data not shown). In 
reality, gene expression dynamics may follow more com- 
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FIG. 3. Effect of r and Pm on cancer cell populations under- 
going prolonged drug treatment, a) Effect of Pm on doubling 
time as a function of r for an initial population of 1000 non- 
mutated cells with gene expression levels above a 95% drug 
threshold. Simulation results using burst model of gene ex- 
pression [15] shown in the inset, b) Heat maps corresponding 
to (a) show probability of remission after 10 generations (100 
realizations). Pm, t and a 2 are scaled by to- 



plex kinetics than that of a simple mean- reverting pro- 
cess due, for example, to multist ability and noise-driven 
switching [I6j [17] . Our simulation results demonstrate 
that such complexity is not required for gene expression 
noise to have an significant impact on population dynam- 
ics under prolonged stress. 

We also determined how the probability of remission 
depends on the mutation rate, the initial number of can- 
cer cells with above-threshold expression, and the time- 
scale of gene expression noise. In these simulations, the 
cancer is in remission if no cells have above threshold gene 
expression and have not acquired a mutation conferring 
permanent immunity within 10 generation times. As ex- 
pected [Fig. ^$jp], the probability of remission is greatly 
decreased when the number of initial surviving cancer 
cells or the mutation rate is increased. Also, when r 
is very short, remission is virtually guaranteed. How- 
ever, the probability that a drug-resistant cell population 
will emerge can be quite substantial within the experi- 
mentally observed range of r. Even with a relative low 
mutation rate (Pm = 0.01) and 10 surviving cells, the 
probability of remission is only 42% when r = 4.0. 

In summary, we have analyzed the effect of gene ex- 
pression noise on the reproductive fitness of isogenic cell 
populations under stress as a first-passage time problem. 
By explicitly incorporating the 'epigenetic memory' of 
this noise (i.e., the fluctuation relaxation time), we have 
generalized previous theoretical work that explained the 
acute effects of noise amplitude but did not incorporate 
gene expression dynamics [5H7]. This generalization is 
important for two reasons. First, it has allowed us to 
demonstrate using a minimal model that gene expression 
noise with biologically realistic time scales has a signif- 
icant effect on reproductive fitness under stress and is 
a critical determinant of population survival. Second, 
it enables theoretical and computational investigations 
of experimentally observed phenomena associated with 
prolonged stress exposure, including reversible shifts in 
gene expression distributions [7], and drug resistance. In 
this context, we have demonstrated that the time-scale of 
epigenetic memory required to develop a drug-resistant 
cell population independently of mutations is compara- 
ble to that measured for certain genes in human cancer 
cells [8]. Correspondingly, long-term population survival 
may not require specialized memory-conferring mecha- 
nisms. It might, for example, be acheived without a 
significant fitness cost through bursty gene expression. 
An important next step is to confirm our findings us- 
ing more realistic models of gene expression incorporat- 
ing additional stochastic effects, such as partitioning er- 
rors [18] , and correspondingly, to employ various analyt- 
ical and numerical methods that may permit solution in 
these more complex cases (e.g. [I9j[20]). We anticipate 
that future analysis of such models will provide a deeper 
understanding of epigenetic interactions between genes, 
drugs and population dynamics. 
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